Figure

This Rmarkdown contain code to obtain all plot necessary to produce our bioinformatical related figures. Final figure anotation and layout were then produce with inkscape.

FIGURE 2

Panel A

# Highlight specific HTO  6x4
Idents(T.Seurat) <- "MULTI_ID"
cellM <- WhichCells(T.Seurat, idents = c("Spleen-M","Thymus-M") )
plotM <- DimPlot(T.Seurat,cells.highlight = cellM,sizes.highlight = 0.5,cols = "grey80", cols.highlight = "coral2",pt.size = 0.5)+ NoLegend()+NoAxes()

cellctrl <- WhichCells(T.Seurat, idents = c("Spleen-ctrl","Thymus-ctrl") )
plotC <- DimPlot(T.Seurat,cells.highlight = cellctrl,sizes.highlight = 0.5,cols = "grey80", cols.highlight="chartreuse4",pt.size = 0.5)+ NoLegend()+NoAxes()

cellMP <- WhichCells(T.Seurat, idents = c("Spleen-MP","Thymus-MP") )
plotMP <- DimPlot(T.Seurat,cells.highlight = cellMP,sizes.highlight = 0.5,cols = "grey80", cols.highlight="#CC0099",pt.size = 0.5)+ NoLegend()+NoAxes()

cellP <- WhichCells(T.Seurat, idents = c("Spleen-P","Thymus-P") )
plotP <- DimPlot(T.Seurat,cells.highlight = cellP,sizes.highlight = 0.5,cols = "grey80", cols.highlight="deepskyblue3",pt.size = 0.5)+ NoLegend()+NoAxes()

#thymus vs rate figure
plotorgan <- DimPlot(T.Seurat, group.by = "tissue",cols=c("darkgoldenrod1","darkred"),pt.size = 0.5)+ NoLegend()+NoAxes()

panA <- grid.arrange(plotM, plotP,plotMP,plotC,plotorgan, nrow = 3)

Panel B

Idents(T.Seurat) <- "integrated_snn_res.1.8"
umapcols <- c("#99A6E7","#E3B1D6","#F78754","#F94A17","#CC79B7","#6477D8","#919296","#CF063A","#FEA981","#55AC7C","#5F4090","#D13A7C","#1332D0","#EC0067","#FB710D","#09B329","#B23ACE","#7A1E45","#FB920D","#FF35AC","#FBAF0D","#ECB90C","#FFDA0D","#E22611")


DimPlot(T.Seurat, reduction = "umap", group.by = "integrated_snn_res.1.8", pt.size = 1,cols = umapcols)+ NoLegend()

panB <- DimPlot(T.Seurat, reduction = "umap", group.by = "integrated_snn_res.1.8", pt.size = 1,cols = umapcols)+ NoLegend()

Panel C

#Dotplot thymus
Idents(T.Seurat.thymus) <- "manualclusters"
T.Seurat.thymus@active.ident <- factor(T.Seurat.thymus@active.ident,levels=c("22","21","20,18,14","8,2,3,23","7","6","10,16","17"))
levels(T.Seurat.thymus)

## [1] "22"       "21"       "20,18,14" "8,2,3,23" "7"        "6"       
## [7] "10,16"    "17"

dotplotthym <- DotPlot(T.Seurat.thymus, dot.scale = 14,features = c("percent.mito","Bmf","Trp53inp1","Tox2","Cd5","Cd69","Cd27","Rag1","Rag2","Cd4","Cd8a","Cd8b1","Mki67","Cdk1","Ptcra","Il2ra","Cd34")) + scale_colour_gradient2(low = "steelblue", mid = "white", high = "red") + coord_flip()+ theme(
   legend.text = element_text(size = 20),
  legend.title = element_text(size = 25),
  axis.text.x = element_text(size = 20,angle = 90),
  axis.text.y = element_text(size = 20),
  axis.title = element_text(size = 25),
  panel.background = element_rect(fill = "grey65",
                                size = 0.5, linetype = "solid"),
  panel.grid.major = element_line(size = 0.5, linetype = 'solid',
                                colour = "grey55"))
#dot scale was 8 but change for rmarkdown output

# Dotplot spleen
Idents(T.Seurat.spleen) <- "manualclusters"

T.Seurat.spleen@active.ident <- factor(T.Seurat.spleen@active.ident,levels=c("10,16","1","4","0,5","12","15","9","17","13","11","19"))
dotplotspleen <- DotPlot(T.Seurat.spleen, dot.scale = 14, features= c("Foxp3","Aes","Anxa1","Gzma","Ccl5","Cxcr3","Sell","S1pr1","Ccr7","Trdc","Tcrg-C4","Tcrg-C2","Tcrg-C1","Trbc2","Trbc1","Trac","Cd8b1","Cd4")) + scale_colour_gradient2(low = "steelblue", mid = "white", high = "red") + theme(
  legend.text = element_text(size = 20),
  legend.title = element_text(size = 25),
  axis.text.x = element_text(size = 20,angle = 90),
  axis.text.y = element_text(size = 20),
  axis.title = element_text(size = 25),
  panel.background = element_rect(fill = "grey65",
                                size = 0.5, linetype = "solid"),
  panel.grid.major = element_line(size = 0.5, linetype = 'solid',
                                colour = "grey55")) + coord_flip()
#dot scale was 8 but change for rmarkdown output

#Combine dotplot
panC <- grid.arrange(dotplotthym,dotplotspleen, ncol = 2)

Panel D

#check par
par(oma = c(0, 0, 0, 0))
par(mar = c(0,1,0,1))
par(mfrow=c(1,2))
#check genotype proportion in each spleen clusters

df <- as.data.frame(as.data.frame.matrix(t(prop.table(table(T.Seurat.spleen@meta.data$HTO,T.Seurat.spleen@meta.data$manualclusters),1)*100)))
df$cluster = rownames(df)
df2<-as.data.frame(t(cbind(rep(60,11),rep(0,11),df)[,1:6]))
rownames(df2[1:2,]) <- c("60","0")

#order data frame
df2 <- df2[c("10,16","19","11","13","17","9","15","12","0,5","4","1")]

radarspleen <- radarchart(df2, cglcol="grey", cglty=1 ,cglwd=0.8, vlcex=0.8, pcol=c("#FF99FF","coral1","cyan3","chartreuse3") , plwd=3, plty=1 ,caxislabels=paste(seq(from = 0,to = 60,by = 15),"%"), axislabcol = "grey40", axistype = 0)

#check genotype proportion in each thymic cluster
df <-  as.data.frame(as.data.frame.matrix(t(prop.table(table(T.Seurat.thymus@meta.data$HTO,T.Seurat.thymus@meta.data$manualclusters),1)*100)))
df$cluster = rownames(df)
df2<-as.data.frame(t(cbind(rep(65,8),rep(0,8),df)[,1:6]))
rownames(df2[1:2,]) <- c("65","0")

#order
df2 <- df2[c("22","17","10,16","6","7","8,2,3,23","20,18,14","21")]


radarthymus <- radarchart(df2, cglcol="grey", cglty=1,caxislabels=paste(seq(0,70,17.5),"%"), axistype = 0,axislabcol="grey40", cglwd=0.8, vlcex=0.8, pcol=c("#FF99FF","coral1","cyan3","chartreuse3") , plwd=3, plty=1 ) 

Figure 4

Panel A

Idents(T.Seurat.spleen) <- "manualclusters"
T.Seurat.spleenbar <- subset(T.Seurat.spleen,  idents = c("0,5","9","1","11","19")) #to only keep cluster needed for the plot
data3 <- data.frame(prop.table(t(prop.table(table(T.Seurat.spleenbar@meta.data$HTO,T.Seurat.spleenbar@meta.data$manualclusters),1)),1)*100)

#order cluster level
data3$Var1 <- factor(data3$Var1, levels = c("0,5","9","1","11","19"))

cellnumber <- data.frame(colSums(table(T.Seurat.spleen@meta.data$HTO,T.Seurat.spleen@meta.data$manualclusters)) )
cellnumber$cluster <- rownames(cellnumber)
row_order <-c("0,5","9","1","11","19")
cellnumber <- cellnumber[row_order,]

Cluster <- data3$Var1
HTO <- data3$Var2
Percentage <- data3$Freq
text <- cellnumber$colSums.table.T.Seurat.spleen.meta.data.HTO..T.Seurat.spleen.meta.data.manualclusters..
# Stacked bar plot
cols <- c("Myc- PTEN- spleen" = "#FF99FF", "MYC- spleen" = "coral1", "PTEN- spleen" = "cyan3", "WT spleen" = "chartreuse3")

ggplot(data3, aes(fill=HTO, y=Percentage, x=Cluster)) + 
    geom_bar(position="stack", stat="identity", width=0.7)+
   xlab("Splenic Clusters")+ylab("Percentage")+ theme(legend.position="bottom")+scale_fill_manual(values =cols, labels=c("Myc Pten","Myc","Pten","WT")) +theme_light()+geom_hline(yintercept=c(25,50,75), linetype="dashed", color = "grey60")+ annotate("text", x = c(1,2,3,4,5), y=103, label = c(text))

Panel B

# On cluster 11 - cd8 memory
C11Ctrl <-rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-ctrl" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "11",])
mC11Ctrl <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C11Ctrl])
C11Pt <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-P" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "11",])
mC11Pt <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C11Pt])
C11MP <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-MP" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "11",])
mC11MP <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C11MP])
C11M <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-M" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "11",])
mC11M <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C11M])
#on cluster 19 - cd8 eff term
C19Ctrl <-rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-ctrl" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "19",])
mC19Ctrl <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C19Ctrl])
C19Pt <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-P" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "19",])
mC19Pt <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C19Pt])
C19MP <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-MP" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "19",])
mC19MP <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C19MP])
C19M <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-M" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "19",])
mC19M <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C19M])
#on cluster 1
C1Ctrl <-rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-ctrl" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "1",])
mC1Ctrl <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C1Ctrl])
C1Pt <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-P" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "1",])
mC1Pt <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C1Pt])
C1MP <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-MP" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "1",])
mC1MP <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C1MP])
C1M <- rownames(T.Seurat.spleen@meta.data[T.Seurat.spleen@meta.data$MULTI_ID == "Spleen-M" & T.Seurat.spleen@meta.data$integrated_snn_res.1.8 == "1",])
mC1M <- mean(T.Seurat.spleen@assays$RNA@data["eYFP",C1M])


l1 <- c("1","1","1","1","11","11","11","11","19","19","19","19")
l2 <- c("Ctrl","Pten","Myc","MycPten","Ctrl","Pten","Myc","MycPten","Ctrl","Pten","Myc","MycPten")
l3 <- c(mC1Ctrl,mC1Pt,mC1M,mC1MP,mC11Ctrl,mC11Pt,mC11M,mC11MP,mC19Ctrl,mC19Pt,mC19M,mC19MP)

tableauYFP <- data.frame(Cluster = l1, Genotypes =l2) 
tableauYFP <- cbind(tableauYFP,Values = l3)


tableauYFP$Cluster <- factor(tableauYFP$Cluster,levels = c("1","11","19"))
tableauYFP$Genotypes <- factor(tableauYFP$Genotypes,levels = c("MycPten","Myc","Pten","Ctrl"))

ggplot(tableauYFP, aes(x = Cluster, Genotypes)) +
        geom_tile(aes(fill = Values)) +
        scale_fill_gradient2( mid='yellow', high='red',limits=c(0,max(tableauYFP$Values)))+ theme_classic(base_size=20)

Supplementary 2

FeaturePlot(T.Seurat, features = "Myc",cols = c("grey", "light blue","cyan3","cyan4","dodgerblue3","blue","mediumslateblue","purple","orchid3","red","brown","black"),order=T)

Supplementary 3

Panel A

Idents(T.Seurat) <- "integrated_snn_res.1.8"
#DGE between our two CD8 naive clusters
dgecd8_data <- read.table(paste(WORKING_DIR, "/02_Seurat_analysis/01_Script/dgeCD8data.txt",sep=""), sep = "\t")
# bar plot
#add abs value to table
dgecd8_data$abs <- abs(dgecd8_data$avg_logFC) 

dgecd8_data$genename <- rownames(dgecd8_data)
# Select markers for plotting on a Heatmap 
markers.use=subset(dgecd8_data, p_val_adj<1e-50 & abs>0.20)
dfcd8markers <-markers.use[order(markers.use$avg_logFC),]

dfcd8markers$genename <- factor(dfcd8markers$genename, levels = dfcd8markers$genename[order(dfcd8markers$avg_logFC)])
dfcd8markers$logpval <- log10(dfcd8markers$p_val_adj)

ggplot(dfcd8markers, aes(x = dfcd8markers$genename, y = dfcd8markers$avg_logFC, fill = logpval)) +   # Fill column
                              geom_bar(stat = "identity", width = .6) +   # draw the bars
                              ylim(-1.2,1.2)+
                              labs(title="DGE - Cluster 1 vs 4 (Cd8 naive)",y ="Log fold change", x = "Genes differentially expressed") +
                              theme_tufte() +  # Tufte theme from ggfortify
                              theme(plot.title = element_text(hjust = .5),axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1), 
                                    axis.ticks = element_blank()) +
                               scale_fill_gradient2(low='red', mid='orange', high='blue',midpoint = -120, breaks=c(-52,-120,-200),labels=c("-50","-120","-200"))+coord_flip() # Flip axes

Panel B

Idents(T.Seurat) <- "integrated_snn_res.1.8"
#DGE between 0 and 12
dgecd4_data <- read.table(paste(WORKING_DIR, "/02_Seurat_analysis/01_Script/dgeCD4data.txt",sep=""), sep = "\t")
# bar plot
#add abs value to table
dgecd4_data$abs <- abs(dgecd4_data$avg_logFC) 

dgecd4_data$genename <- rownames(dgecd4_data)
# Select markers for plotting on a Heatmap 
markers.use=subset(dgecd4_data,p_val_adj<1e-10 & abs>0.2)
dfcd4markers <-markers.use[order(markers.use$avg_logFC),]

dfcd4markers$genename <- factor(dfcd4markers$genename, levels = dfcd4markers$genename[order(dfcd4markers$avg_logFC)])
dfcd4markers$logpval <- log10(dfcd4markers$p_val_adj)

ggplot(dfcd4markers, aes(x = dfcd4markers$genename, y = dfcd4markers$avg_logFC, fill = logpval)) +   # Fill column
                              geom_bar(stat = "identity", width = .6) +   # draw the bars
                              ylim(-1.2,1.2)+
                              labs(title="DGE - Cluster 0 vs 12 (Cd4 naive)",y ="Log fold change", x = "Genes differentially expressed") +
                              theme_tufte() +  # Tufte theme from ggfortify
                              theme(plot.title = element_text(hjust = .5),axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1), 
                                    axis.ticks = element_blank()) +
                               scale_fill_gradient2(low='red', mid='orange', high='blue',midpoint = -40, breaks=c(-12,-40,-57),labels=c("-10","-40","-60"))+ coord_flip() # Flip axes

Panel C

Idents(T.Seurat) <- "integrated_snn_res.1.8"

# get all gene name express in our cells as background
background <- T.Seurat@assays$RNA@meta.features
backgroundrow <- rownames(background)

genecomprow <- read.table(paste(WORKING_DIR, "/02_Seurat_analysis/01_Script/cluster1v4.txt",sep=""), sep = "\t")
genecomprow$x = as.character(genecomprow$x)
genecomprow <- genecomprow[,1]

CPenrich <- enrichGO(gene= genecomprow, OrgDb = 'org.Mm.eg.db', ont="BP",keyType = "SYMBOL",universe = backgroundrow) # org.Mm.eg.db genome mouse
#head (CPenrich)

dotplot(CPenrich, showCategory=15,color = "p.adjust",x="count") #+ coord_flip()+theme(axis.text.x = element_text(angle = 90, hjust = 1))

Panel D

genecomprow <- read.table(paste(WORKING_DIR, "/02_Seurat_analysis/01_Script/cluster0vs12.txt",sep=""), sep = "\t")
genecomprow$x = as.character(genecomprow$x)
genecomprow <- genecomprow[,1]
####


CPenrich <- enrichGO(gene= genecomprow, OrgDb = 'org.Mm.eg.db', ont="BP",keyType = "SYMBOL",universe = backgroundrow) # org.Mm.eg.db genome mouse

dotplot(CPenrich, showCategory=15,color = "p.adjust",x="count")+ scale_y_discrete(labels=function(x)str_wrap(x, width=40))

Supplementary 4

Panel A

#DOT plot eyFP and TGD on CD8 effector and memory
Idents(T.Seurat.spleen) <- "integrated_snn_res.1.8"
CD8sub <- subset(T.Seurat.spleen, idents = c("11","19","13"))
Idents(CD8sub) <- "MULTI_ID"
CD8sub@active.ident <- factor(CD8sub@active.ident,levels=c("Spleen-M","Spleen-MP","Spleen-ctrl","Spleen-P"))
#levels(CD8sub)

DotPlot(CD8sub, dot.scale = 8,features = c("Tcrg-C1","Trdc","Trbc2","Trac","eYFP") ) + scale_colour_gradient2(low = "steelblue", mid = "white", high = "red")+ ggtitle("CD8 memory and effector clusters")